/*

 ARCHIVO: "complex.h"

 Declaración de la clase complejo.

 Primera versión: 23-3-2012

*/


#include <cmath>
#include <iostream>

#include "globals.h"
#include "complex.h"


/* ------------------------------- Definición de la clase complejo - Begin ------------------------------- */


//Constructor en base a dos doubles, pudiendo ser interpreatdos en forma cartesiana o polar
//según el valor de bpolar.
complejo::complejo(double da, double db, bool bpolar) : re_(0.0), im_(0.0)
{
    if (bpolar)
        re_ = da*std::cos(db), im_ = da*std::sin(db);
    else
        re_ = da, im_ = db;
}

//Constructor copia.
complejo::complejo(const complejo &cmr) : re_(cmr.re_), im_(cmr.im_)
{
}

//Destructor.
complejo::~complejo()
{
}

//Parte real del complejo.
double complejo::real() const
{
    return re_;
}

//Parte imaginaria del complejo.
double complejo::imag() const
{
    return im_;
}

//Módulo del complejo.
double complejo::abs() const
{
    return sqrt(re_ * re_ + im_ * im_);
}

//Módulo del complejo al cuadrado, igual al complejo por su conjugado.
double complejo::abs_cuad() const
{
    return re_ * re_ + im_ * im_;
}

//Argumento del complejo.
double complejo::arg() const
{
    if (re_ > 0.0)
        return std::atan(im_/re_);
    else if (re_ < 0.0)
    {
        if (im_ > 0.0)
            return std::atan(im_/re_) + M_PI;
        else if (im_ < 0.0)
            return std::atan(im_/re_) - M_PI;
        else
            return M_PI;
    }
    else
    {
        if (im_ > 0.0)
            return M_PI_2;
        else if (im_ < 0.0)
            return -M_PI_2;
        else
            return 0.0; // Arbitrario, no está definido.
    }
}

//Conjuga la instancia del complejo.
complejo& complejo::conjugar()
{
    im_ *= -1.0;
    return *this;
}

//Devuelve otro complejo que es el conjugado de la instancia de complejo.
complejo complejo::conjugado() const
{
    complejo cmr(*this);
    cmr.im_ *= -1.0;
    return cmr;
}


/* Operadores - Begin */

//Operador asignación (complejo).
complejo& complejo::operator=(const complejo &cmr)
{
    re_ = cmr.re_;
    im_ = cmr.im_;
    return *this;
}

//Operador asignación (double).
complejo& complejo::operator=(double cmr)
{
    re_ = cmr;
    im_ = 0.0;
    return *this;
}

//Operador suma con asignación.
complejo& complejo::operator+=(const complejo &cmrightop)
{
    re_ += cmrightop.re_, im_ += cmrightop.im_;

    return *this;
}

//Operador resta con asignación.
complejo& complejo::operator-=(const complejo &cmrightop)
{
    re_ -= cmrightop.re_, im_ -= cmrightop.im_;

    return *this;
}

//Operador producto con asignación.
complejo& complejo::operator*=(const complejo &cmrightop)
{
    double re = re_ * cmrightop.re_ - im_ * cmrightop.im_;
    double im = re_ * cmrightop.im_ + im_ * cmrightop.re_;

    re_ = re, im_ = im;

    return *this;
}

//Operador división con asignación.
complejo& complejo::operator/=(const complejo &cmrightop)
{
    double modquad = cmrightop.abs_cuad();

    complejo rconj = cmrightop.conjugado();

    double re = (re_ * rconj.re_ - im_ * rconj.im_)/modquad;
    double im = (re_ * rconj.im_ + im_ * rconj.re_)/modquad;

    re_ = re, im_ = im;

    return *this;
}

//Operador menos unario.
complejo operator-(const complejo &cmr)
{
    return complejo(-cmr.real(), -cmr.imag());
}

//Operador suma.
complejo operator+(const complejo& cmleftop, const complejo& cmrightop)
{
    complejo aux = cmleftop;

    return (aux += cmrightop);
}

//Operador resta.
complejo operator-(const complejo &cmleftop, const complejo &cmrightop)
{
    complejo aux = cmleftop;

    return (aux -= cmrightop);
}

//Operador producto.
complejo operator*(const complejo &cmleftop, const complejo &cmrightop)
{
    complejo aux = cmleftop;

    return (aux *= cmrightop);
}

//Operador división.
complejo operator/(const complejo &cmleftop, const complejo &cmrightop)
{
    complejo aux = cmleftop;

    return (aux /= cmrightop);
}

//Operador igualdad.
bool operator==(const complejo &cmleftop, const complejo &cmrightop)
{
    return ((cmleftop.real() == cmrightop.real()) && (cmleftop.imag() == cmrightop.imag()));
}

//Operador distinto.
bool operator!=(const complejo &cmleftop, const complejo &cmrightop)
{
    return ((cmleftop.real() != cmrightop.real()) || (cmleftop.imag() != cmrightop.imag()));
}

//Operador de salida por stream.
std::ostream& operator<<(std::ostream &os, const complejo &cmr)
{
    return (os << "(" << cmr.real() << "," << cmr.imag() << ")");
}

//Operador de entrada desde stream.
std::istream& operator>>(std::istream &is, complejo &cmr)
{
    double re = 0.0;
    double im = 0.0;
    char ch = '\0';
    bool okay = false;

    if(is >> ch)
    {
        if ((ch == '(') && (is >> re)  && (is >> ch) &&
                (ch == ',') && (is >> im) && (is >> ch) && (ch == ')'))
            okay = true;
        else
        {
            is.putback(ch);

            if (is >> re)
                okay = true;
        }
    }

    if (okay)
        cmr = complejo(re, im);
    else
        is.clear(std::ios::badbit);

    return is;
}

/* Operadores - End */

/* ------------------------------- Definición de la clase complejo - End ------------------------------- */




/* Funciones asociadas - Begin */

//Funcción que calcula la exponencial de un complejo.
complejo exp(const complejo &cmarg)
{
    return complejo(std::exp(cmarg.real()), cmarg.imag(), true);
}

//Función que calcula la potencia de un número double elevado a un complejo.
complejo pow(double dr, complejo cmexp)
{
    return complejo(std::pow(dr, cmexp.real()), std::log(dr)*cmexp.imag(), true);
}

//Funciones que calculan la potencia de un número complejo.
complejo pow(const complejo &cmr, int iexp)
{
    return complejo(std::pow(cmr.abs(), iexp), cmr.arg()*iexp, true);
}

complejo pow(const complejo &cmr, double dexp)
{
    return complejo(std::pow(cmr.abs(), dexp), cmr.arg()*dexp, true);
}

complejo pow(const complejo &cmr, complejo cmexp)
{
    return pow(cmr.abs(), cmexp)*exp(cmr.arg()*cmexp*complejo(0, 1));
}

//Función que calcula la dft de un vector de complejos.
//Usada para verificar la correctitud de las operaciones.
void dft(complejo *in, complejo *out, int n)
{
    for (int i = 0 ; i < n ; ++i)
    {
        out[i] = 0;

        for (int j = 0 ; j < n ; ++j)
        {
            complejo arg(0, -2*M_PI*i*j/double(n));
            out[i] = out[i] + in[j] * exp(arg);
        }
    }
}

// Función de testeo de la aritmética de complejos
// usando la DFT.
// El resultado se compara luego con lo obtenido
// con el comando fft de MATLAB.

#define LEN 20

void test_complejo()
{
    complejo a[LEN], b[LEN];

    for (int i = 0; i < LEN; i++)
        a[i] = i;

    dft(a, b, LEN);

    for (int i = 0; i < LEN; i++)
        std::cout << b[i];
}

/* Funciones asociadas - End */


